/* Integration points (quadratures)
   Jukka Toivanen, Jan Stebel
*/

#include "linbubblepoint.h"



template<>
void LinBubblePoint<3>::setSF()
{
  SF(0,0) = 1-ksi[0]-ksi[1]-ksi[2];
  SF(1,0) = ksi[0];
  SF(2,0) = ksi[1];
  SF(3,0) = ksi[2];
  SF(4,0) = ksi[0]*ksi[1]*ksi[2]*(1-ksi[0]-ksi[1]-ksi[2]);
}

template<>
void LinBubblePoint<3>::setGSF()
{
  GSF(0,0) = -1;
  GSF(0,1) = 1;
  GSF(0,2) = 0;
  GSF(0,3) = 0;
  GSF(0,4) = ksi[1]*ksi[2]*(1-2*ksi[0]-ksi[1]-ksi[2]);
  GSF(1,0) = -1;
  GSF(1,1) = 0;
  GSF(1,2) = 1;
  GSF(1,3) = 0;
  GSF(1,4) = ksi[0]*ksi[2]*(1-ksi[0]-2*ksi[1]-ksi[2]);
  GSF(2,0) = -1;
  GSF(2,1) = 0;
  GSF(2,2) = 0;
  GSF(2,3) = 1;
  GSF(2,4) = ksi[0]*ksi[1]*(1-ksi[0]-ksi[1]-2*ksi[2]);
}

template<>
inline int LinBubblePoint<3>::mapDof(int elem, int n)
{
  return (n==4)?(imesh->getNnodes()+elem):(imesh->getVolumeNode(elem, n));
}

template<>
inline int LinBubblePoint<3>::mapDof(int elem, int dof, int iParent)
{
	throw 1;
}







template<>
void LinBubblePoint<2>::setSF()
{
  SF(0,0) = 1-ksi[0]-ksi[1];
  SF(1,0) = ksi[0];
  SF(2,0) = ksi[1];
  SF(3,0) = ksi[0]*ksi[1]*(1-ksi[0]-ksi[1]);
}

template<>
void LinBubblePoint<2>::setGSF()
{
  GSF(0,0) = -1;
  GSF(0,1) = 1;
  GSF(0,2) = 0;
  GSF(0,3) = ksi[1]*(1-2*ksi[0]-ksi[1]);
  GSF(1,0) = -1;
  GSF(1,1) = 0;
  GSF(1,2) = 1;
  GSF(1,3) = ksi[0]*(1-ksi[0]-2*ksi[1]);
}

template<>
inline int LinBubblePoint<2>::mapDof(int elem, int n)
{
  return (n==3)?(imesh->getNlnodes()+elem):(imesh->getFaceNode(elem, n));
}

template<>
inline int LinBubblePoint<2>::mapDof(int elem, int dof, int iParent)
{
	return (dof==3)?-1:imesh->getVolumeNode(imesh->getFaceParentNo(elem,iParent),imesh->getFaceParentNodeNo(elem,iParent,dof));
}



template<>
void LinBubblePoint<1>::setSF()
{
  SF(0,0) = 1-ksi[0];
  SF(1,0) = ksi[0];
  SF(2,0) = ksi[0]*(1-ksi[0]);
}

template<>
void LinBubblePoint<1>::setGSF()
{
  GSF(0,0) = -1;
  GSF(0,1) = 1;
  GSF(0,2) = 1-2*ksi[0];
}

template<>
inline int LinBubblePoint<1>::mapDof(int elem, int n)
{
  return (n==2)?(imesh->getNlnodes()+elem):(imesh->getEdgeNode(elem, n));
}

template<>
inline int LinBubblePoint<1>::mapDof(int elem, int dof, int iParent)
{
	return (dof==2)?-1:imesh->getFaceNode(imesh->getEdgeParentNo(elem,iParent),imesh->getEdgeParentNodeNo(elem,iParent,dof));
}






template class LinBubblePoint<1>;
template class LinBubblePoint<2>;
template class LinBubblePoint<3>;


